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Summary 

A three-dimensional cubic cavity flow has been analyzed for diatomic 
gases by using the Boltzmann equation with the Bhatnagar-Gross-Rrook 
(B-G-K) model. The method of discrete ordinate was applied, and the 
diffuse reflection boundary condition was assumed. The results, which show 
a consistent trend toward the Navier-Stokes solution as the Knudson 
number is reduced, give us confidence to apply the method to a three- 
dimensional geometry for practical predictions of rarefied-flow 
characteristics. The CPU time and the main memory required for a three- 
dimensional geometry using this method seem reasonable. 

Symbols 

a defined in equations (5) and (6) 

energy level (equal to the roots of the Laguerre polynomial 
of degree n) 

F Maxwellian distribution function 

f distribution function 

imax maximum i index 

jmax maximum j index 

Kn Knudson number 

k Boltzmann constant 

kmax maximum k index 
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weighting coefficient of the odd equally spaced quadrature 
length 

mass of a molecule 
number particle density 
gas constant 

Gauss-Laguerre weighting coefficient 
temperature 

macroscopic velocity (flow velocity) 
molecular velocity 
velocity of moving surface 
vertical distance 

rotational relaxational parameter 

mean free path 

viscosity 

collision frequency 
nondimensional quantity 

elastic 

equilibrium 

index of physical space 

inelastic 

total 

translational 

x-component 

y-component 

z-component 

index of molecular velocity space 
discrete energy level 
reference condition 


Introduction 

With space vehicles orbiting the world at hypersonic speed, numerical 
methods to cover all flow regimes become more important. Advances in 
computer hardware in recent years allow methods based on the kinetic 
theory to be used for practical applications. Many two-dimensional 
problems have been solved by the kinetic approach. A limited number 
of three-dimensional problems have been solved recently by the direct- 
simulation Monte Carlo (DSMC) method . 8 ' 10 However, there are no three- 
dimensional numerical solutions based on the discrete ordinate method 11 
applied to a Boltzmann equation. 

A cubic cavity flow problem has been solved for diatomic gases using 
the Boltzmann equation with the Bhatnagar-Gross-Krook (B-G-K) 
model . 12 ’ 13 Because neither theoretical nor experimental data are available 
for cubic cavity flow in free molecule flow, slip How, and transition flow 
regimes, efforts have been made to obtain a continuum flow solution as 
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closed as possible so that the results from this study can be compared with 
the results of the Navier-Stokes solution. 14 

Governing Equation 

The geometry of the problem considered is shown in figure 1. A cubic 
cavity with the same dimension on all sides is filled with diatomic gas (such 
as air), and the wall temperature is assumed to be the same on all walls. 
The top surface (j = jmax surface) is moving at a constant speed, and 
steady-state conditions are reached and analyzed numerically. 

The distribution function at the energy level E ff , f , is governed by the 
Boltzmann equation with the Bhatnagar-Gross-Krook-Morse model for the 
collision terms: 
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The elastic collision frequency t / e | 15,16 is given by 
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whereas the 

inelastic collision frequency i/. . is related to v , 

in the 

formula of 
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where the rotational relaxational parameter Z R must be obtained by 
experiment. The viscosity-temperature relation of the form 
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( 7 ) 


Ml 

is assumed and the value of s can be found in reference 17; for the current 
study, a and s are equal to 0.4 and 0.756, respectively, for air. The 
subscript 1 indicates the reference condition. The reference condition used 
in this study is the equilibrium flow in the cubic cavity when the top surface 
is stationary. 

The reference viscosity is related to the reference mean free path Aj 
by the relation 



= _Lmn 1 A 1 (2jrRT 1 ) 1 / 2 
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Combining equations (4) to (8) gives 
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All macroscopic quantities, which will be listed in the nondimensional 
form, can be calculated in terms of f . 

1 (T 

Nondimensionalization 

A characteristic velocity used for nondimensionalization is defined as 



V x = (2RTJ) 1 / 2 (10) 

The definitions of nondimensional variables are defined as follows: 
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Equation (1) becomes 
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Discretization of Distribution Function and Numerical Scheme 

The distribution function at energy level a, f^., was discretized with 
equal spacing in the physical space, the arguments of the odd equally spaced 
quadrature (see appendix and ref. 18) was the spacing for the velocity space, 
and the roots of the Laguerre polynomial of degree 4 was used for the 
spacing for the internal energy space. The value of the distribution function 
at these discretized points was calculated by solving the Boltzmann equation 
using the finite difference method in the physical space. 

The difference scheme is the first-order upwind difference based on the 
molecular velocity; for example, 
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The distribution functions were then integrated for all macroscopic 
properties by using the odd equally spaced quadrature, e.g., the particle 
density at energy level a and at the physical point (i, j, k) is 
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where k a , k^, and k^ are the weighting coefficients of the odd equally 
spaced quadrature for the velocity components v a , v^, and v^, respec- 
tively. These properties were summed over the energy space, and then 
saved and used for the next iteration. 


Boundary Conditions 

For a constantly moving surface, perfectly diffuse reflection is assumed to 
specify the interaction of the molecules with the surface of the moving plate. 
Molecules which strike the moving surface (j = jmax surface) are emitted 
with a Maxwellian velocity distribution characterized by the plate temper- 
ature T w and the plate velocity w; that is, 
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where the density of molecules diffusing from the plate n w is not known 
a priori and must be found by applying the condition of zero mass flux 
normal to the plate at the surface, i.e. , 
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where f is the distribution function for v„ > 0 and f is the distri- 
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bution function for v y < 0. In the discrete ordinate form, this gives the 
relation at (i, j, k) 
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Similarly, a zero-mass-flux boundary condition was applied to all other 
surfaces, and a Maxwellian velocity distribution function characterized by 
the plate temperature was used for the molecules emitted from the surface. 

Because of the symmetry of the flow inside the cubic cavity, half of the 
cubic cavity was used for the computation. At the symmetry plane (i.e., 
k = kmax plane), the symmetry boundary condition was applied; that is, 
the distribution function of the outgoing molecule is equal to the distri- 
bution function of the incoming molecule at the symmetry plane. 


Computational Procedures 

For the present study, the initial condition assumed was that the flow 
everywhere was characterized by n = 1 and T = 1 . As mentioned 
previously, the j = jmax surface was selected to be the moving surface. 
The computation began from the j = 2 surface with a known value for 
at the j = 1 surface. On each j = constant surface, the computation 
was performed along the i direction first and then marching toward the 
symmetry plane (k = kmax) starting with f^ + The symbol f (7 + + + is 
defined as the subdistribution function at energy level <r for v x > 0, 
v y > 0, and v z > 0. At the symmetry plane, the symmetry boundary con- 
dition was applied (i.e., f = f + f ). Knowing f^ at the sym- 
metry plane, the computation was marched backward from k = kmax - 1 
toward k = 1. Similarly, the computation continued for f "* ^ and f^ + . 

After the computation reached the j = jmax surface, the particle 
density was found at this surface by applying the no-flux boundary 
condition and then the computation marched downward from j = jmax - 1 
to j = 1 surface. The order of computation was f^ + + , 

t a , f^. , and f^ . The particle density on the surfaces could then 

be found by using the no-flux boundary condition and the iteration 


8 


continued. The fourth-order Gauss-Laguerre quadrature was used for the 
energy space, and four sets of the third-order odd equally spaced quadrature 
were used for the molecular velocity space with the equal spacing of 0.15. 
(See the appendix.) The convergence was assumed when the maximum 
particle density increment was less than 10” 4 . 

Results 

As mentioned previously, efforts have been made to calculate the flow 
with as small a Knudson number as possible so that the results could be 
compared with the solution of the Navier-Stokes equation. The results of 
this study are all based on w = 0.1. Comparisons between the Boltzmann 
solution and the Navier-Stokes solution 14 are given in Figures 2 to 4. 

Figure 2 shows the velocity profiles along the centerline for free 
molecule flow (Knudson number (Kn) = 100), slip flow (Kn = 0.1), and 
nearly continuum flow (Kn = 0.03). It clearly shows that the flow slips on 
the moving surface for Kn = 100 and Kn = 0.1; that is, the flow velocity 
on the surface is less than the speed of the moving wall. It also shows that 
the no-slip boundary condition of the Navier-Stokes equation was almost 
recovered for Kn = 0.03. 

The comparison of the velocity vector plot on the symmetry plane 
between the Navier-Stokes solution and the Boltzmann solution is given in 
figure 3. The general shapes of the primary vortex are similar except that 
the center of the vortex for the Navier-Stokes solution is slightly farther 
downstream than that for the Boltzmann solution. 

Figure 4 shows the surface static pressure looking from the center of the 
cubic cavity toward the upper corner for both the Navier-Stokes solution 
and the Boltzmann solution. The three-dimensional effect of the static 
pressure distribution is clearly shown qualitatively. The Boltzmann 
solution, even for the case of Kn = 0.03, is still not close to the continuum 
solution obtained by the Navier-Stokes equation; however, the consistent 
trend toward the Navier-Stokes solution is encouraging. 

The symmetry plane velocity vector plots for two different Knudson 
numbers are shown in figure 5. The shape of the primary vortex is clearly 
shown even for Kn = 100. The center of the vortex moves upward toward 
the center obtained by the Navier-Stokes solution (fig. 3) when the Knudson 
number is reduced. It also clearly shows that the magnitude of the velocity 
vectors increases as expected when the Knudson number is reduced. 

Figure 6 is a plot of the number particle density on the surfaces for 
three Knudson numbers. The distribution patterns for Kn = 100 and 
Kn = 0.03 are completely different. For a free molecular flow (Kn = 100), 
there are not enough molecular collisions to ensure the high number particle 
density at the upper downstream corner as shown for the case of Kn = 0.03. 
The high number particle density on the upper corner of the downstream 
vertical surface (i.e. , near i = imax and j = jmax) and the low number 
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particle density on the upper corner of the upstream vertical surface (i.e., 
near i = 1 and j = jmax) is a reasonable solution. 

The convergence history is shown in figure 7 for Kn = 0.03. Each 
iteration took about 40 min on the Cray 2 supercomputer. Fortunately, it 
required less than 50 iterations to reach a convergent solution. In other 
words, it took less than 30 hr of CPU time for Kn = 0.03. The CPU time 
for Kn = 100 was less than 10 hr. To reduce the main memory 
requirement, the Maxwellian distribution functions and collision frequency 
were calculated repeatedly; as a consequence, extra CPU time was needed 
for each iteration. The grid size used in the study was 31 by 31 by 16 in 
the physical space, 12 by 12 by 12 in the molecular velocity space and 4 
levels of internal energy in the energy space. Considerable CPU time can 
be saved by improving the integration scheme in the molecular velocity 
space. In any case, the discrete ordinate method used in this study can be 
applied to simple three-dimensional geometries without using the parallel 
programming technique. 


Concluding Remarks 

A three-dimensional cubic cavity flow was solved for diatomic gases by 
using the Boltzmann equation with the Bhatnagar-Gross-Krook (B-G-K) 
model. A comparison was made between the Boltzmann solution and the 
Navier-Stokes solution for the velocity profiles along the centerline, the 
primary vortex on the symmetry plane, and the surface static pressure. The 
general trend toward the Navier-Stokes solution as the Knudson number is 
reduced indicates that the solutions are very reasonable and that the 
discrete ordinate method can be used with confidence to a three-dimensional 
geometry for practical predictions of rarefied-flow characteristics. Because 
of the robustness of this numerical scheme, it requires less than 50 iterations 
to obtain a converged solution. The present method has a potential to be 
a practical flow simulation method to cover all flow regimes. 

Appendix — Laguerre and Odd Equally Spaced Quadratures 
Abscissa and Weight Factors for Laguerre Integration 


Je- E f(E)dE = £ R^E,) 
0 i- 1 


where E ; are the abscissas and are the roots of Laguerre polynomials. 
For n = 4, 

Ej = 0.322547689619 
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E 2 = 1.745761101158 
E 3 = 4.536620296921 
E 4 = 9.395070912301 
= 6. 03154104342x10" 1 
R 2 = 3.57418692438 x 10" 1 
R 3 = 3. 88879085150 xlO" 2 
R 4 = 5.39294705561 XlO" 4 

Odd Equally Spaced Quadrature 

1 n 

Jf(x)dx ~ £ V(“i) 
0 i=l 


where the arguments are taken to be l/(2n), 3/(2n),..., 1 - [l/(2n)] and 
the weighting coefficients kj for n = 3 are 

k x = k 3 = 0.375 
k 2 = 0.25 
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w = 0.1 


yfeJ 

*— Centerline 


Figure 2. — Velocity profiles along centerline. 





(a) Navier-Stokes solution (Re = 1 00). (b) Boltzmann solution (Kn = 0.03). 

Figure 3.— Velocity vectors on symmetry plane for Navier-Stokes solution and Boltzmann solution. 
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(a) Navier-Stokes solution (Re = 1 00). (b) Boltzmann solution (Kn = 0.03). 

Figure 4. — Surface static pressure for Navier-Stokes solution and Boltzmann solution. 



(a) Kn = 1 00. (b) Kn = 0.03. 

Figure 5.— Velocity vectors on symmetry plane for two Knudson numbers. 
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(a) Kn = 0.1 . 



(b) Kn = 1 00. (c) Kn = 0.03. 

Figure 6. — Number particle density on surfaces for three Knudson numbers. 
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